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Abstract 

Wc construct finite-dimensional approximations of solution spaces of divergence 
form operators with ^^-coefficients. Our method does not rely on concepts of er- 
godicity or scale-separation, but on the property that the solution space of these 
operators is compactly embedded in H 1 if source terms are in the unit ball of L 2 in- 
stead of the unit ball of H~ x . Approximation spaces are generated by solving elliptic 
PDEs on localized sub-domains with source terms corresponding to approximation 
bases for H 2 . The 7J 1 -error estimates show that C(/i~ d )-dimcnsional spaces with 
basis elements localized to sub-domains of diameter 0(h a In i) (with a € [|, 1)) 
result in an 0(h 2 ~ 2a ) accuracy for elliptic, parabolic and hyperbolic problems. For 
high-contrast media, the accuracy of the method is preserved provided that local- 
ized sub-domains contain buffer zones of width 0(h a In ^) where the contrast of the 
medium remains bounded. The proposed method can naturally be generalized to 
vectorial equations (such as elasto-dynamics) . 

1 Introduction 

Consider the partial differential equation 



- div i^a(x)Vu(x)J = g(x) x £ 0; g £ L 2 (Q), a(x) = {a„ G L°°(ft)} 
u = on <9f2, 



(1.1) 



where O is a bounded subset of M d with a smooth boundary (e.g., C 2 ) and a is symmetric 
and uniformly elliptic on Q. It follows that the eigenvalues of a are uniformly bounded 
from below and above by two strictly positive constants, denoted by A m i n (a) and A max (a). 
Precisely, for all £ E M. d and x£ll, 

A mi n(a)|£| 2 < £ T a(*)£ < A max (a)|e| 2 . (1.2) 
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In this paper, we are interested in the homogenization of (jl.ip (and its parabolic 
and hyperbolic analogues in Sections H] and [5]) , but not in the classical sense, i.e., that 
of asymptotic analysis [9] or that of G or .ff -convergence ([17], |57j 132] ) in which one 
considers a sequence of operators — div(a e V) and seeks to characterize limits of solution. 
We are interested in the homogenization of (jl.ip in the sense of "numerical homogeniza- 
tion," i.e., that of the approximation of the solution space of (jl.ip by a finite-dimensional 
space. 

This approximation is not based on concepts of scale separation and/or of ergodicity 
but on compactness properties, i.e., the fact that the unit ball of the solution space 
is compactly embedded into H^iQ) if source terms (g) are integrable enough. This 
higher integrability condition on g is necessary because if g spans .ff _1 (f2), then the 
solution space of (jl.ip is Hq(£1) (and it is not possible to obtain a finite dimensional 
approximation subspace of Hq(Q) with arbitrary accuracy in i/ 1 -norm). However, if 
g spans the unit ball of L 2 (fl), then the solution space of (jl.ip shrinks to a compact 
subset of Hq(Q.) that can be approximated to an arbitrary accuracy in ff 1 -norm by 
finite-dimensional spaces [10] (observe that if a = Id, then the solution space is a closed 
bounded subset of H 2 n Hq (O), which is known to be compactly embedded into Hq(Q)). 

The identification of localized bases spanning accurate approximation spaces relies 
on a transfer property obtained in |10j . For the sake of completeness, we will give a 
short reminder of that property in Section [2j In Section [3j we will construct localized 
approximation bases with rigorous error estimates (under no further assumptions on 
a than those given above) . In Sub-section 13. 4j we will also address the high-contrast 
scenario in which A max (a) is allowed to be large. In Sections H] and [SJ we will show that 
the approximation spaces obtained by solving localized elliptic PDEs remain accurate for 
parabolic and hyperbolic time-dependent problems. We refer to Section [6] for numerical 
experiments. We refer to Section [B] of the Appendix for further discussion and a proof 
of the strong compactness of the solution space when the range of g is a closed bounded 
subset of H~ U (Q) with v < 1 (this notion of strong compactness constitutes a simple 
but fundamental link between classical homogenization, numerical homogenization and 
reduced order modeling). 

2 A reminder on the flux-norm and the transfer property. 

Recall that the key element in G and H convergence is a notion of "compactness by 
compensation" combined with convergence of fluxes. Here, the notion of compactness is 
combined with a flux- norm introduced in |10j . 

The flux-norm. We will now give a short reminder on the flux-norm and its properties. 

Definition 2.1. For k G (L 2 (Q,)) d , denote by k po t the potential portion of the Weyl- 
Helmholtz decomposition of k. Recall that k po t is the orthogonal projection of k onto 
{V/ : feH^)}m(L 2 (n)) d . 
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Definition 2.2. For ip G H^(U), define 

H^IU-flux := II (a V^)pot || (L2 (n)) d. (2.1) 
We call 1 1 I la-flux the flux- norm of VP. 



The following proposition shows that the flux-norm is equivalent to the energy norm 
if A m i n (a) > and A mui (a) < oo. 

Proposition 2.1. [Proposition 2.1 of [10]] ||.|| a -flux ^ s a norm on Hq(Q). Furthermore, 
for all if) G ^q 1 ^) 

Amm(a)||V^|| (L2(n)) d < ||^||a-flux < A max (a)|| V^|| (i2(a ))d. (2.2) 

Motivations behind the flux-norm: There are three main motivations behind the 
introduction of the flux norm. 

• The flux-norm allows to obtain approximation error estimates independent from 
both the minimum and maximum eigenvalues of a. In fact, the flux- norm of the 
solution of (jl.ip is independent from a altogether since 



u 



-flux = || V A 1 fi'||(L2( n ))d. (2.3) 



• The (-)pot in the a-flux-norm is explained by the fact that in practice, we are inter- 
ested in fluxes (of heat, stress, oil, pollutant) entering or exiting a given domain. 
Furthermore, for a vector field £, Jqq^ • nds = J n div(^ pot )dx, which means that 
the flux entering or exiting is determined by the potential part of the vector field. 

• Classical homogenization is associated with two types of convergence: convergence 
of energies (r-convergence [331115] ) and convergence of fluxes (G or iJ-convergence 
|47t 152"! [5*5] 157] 146]). Similarly, one can define an energy norm and a flux- norm. 

The transfer property. For V, a finite dimensional linear subspace of Hq(£1), we 
define 

(divaW) := span{div(aVi;) : v G V}. (2.4) 
Note that (divaV^) is a finite dimensional subspace of H~ 1 (Q,). 

Theorem 2.1. (Transfer property of the flux norm) [Theorem 2.1 of [ID] ] Let V 

and V be finite- dimensional subspaces of Hq(Q). For f G L 2 (U), let u be the solution 
°f (COP with conductivity a and u' be the solution of M-l\) with conductivity a'. If 
(divaW) = (div a'W); then 

inf Ik - HU-flux = ^ IK - HU'-flux ^ (2 5) 

v ^ v llfflk 2 ^) vev \\g\\L 2 (n) 
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The usefulness of (|2.5|) can be illustrated by considering a! = I so that diva'V = A. 
Then, v! E H 2 and therefore V can be chosen as, e.g., the standard piecewise linear 
FEM space, on a regular triangulation of of resolution h, with nodal basis {4>i}. The 
space V is then defined by its basis {Oi} determined by 

div(aWi) = Afa in n 
6i = on dn. 

Equation (|2.5p shows that the approximation error estimate associated with the space 
V and the problem with arbitrarily rough coefficients is (in a-flux norm) equal to the 
approximation error estimate associated with piecewise linear elements and the space 
H 2 (Q,). More precisely, 

sup inf ^ ~ ^"- flux < C h, (2.7) 

where C does not depend on a. 

We refer to [22], [25] and for recent results on finite element methods for high 
contrast (A max (a)/A m i n (a) >> 1) but non-degenerate (A m i n (a) = 0(1)) media under 
specific assumptions on the morphology of the (high-contrast) inclusions (in [22], the 
mesh has to be adapted to the morphology of the inclusions). Observe that the pro- 
posed method remains accurate if the medium is both of high contrast and degenerate 
(A m j n (a) << 1), without any further limitations on a, at the cost of solving PDEs (|2.6p 
over the whole domain S7. 

Remark 2.1. We refer to [10] for the optimal constant C in (|2.7p . This question of 
optimal approximation with respect to a linear finite dimensional space is related to the 
Kolmogorov n- width [544 114"] , which measures how accurately a given set of functions 
can be approximated by linear spaces of dimension n in a given norm. A surprising 
result of the theory of n-widths is the non-uniqueness of the space realizing the optimal 
approximation [55]. Observe also that, as another consequence of the transfer property 
(|2.5|) . a h k+1 rate of convergence can be achieved in (|2.7p by replacing fa with higher- 
order basis functions in (|2.6p . and ||g||£2 with ||<7||^fc in (|2.7p . Similarly an exponential 
rate of convergence can be achieved if the source terms g are analytic. This is the 
reason behind the near exponential rate of convergence observed in [6] for harmonic 
functions (i.e., with zero source terms, and particular "buffer" solutions computed near 
the boundary) and bounded (non high) contrast media. 

3 Localization of the transfer property. 

The elliptic PDEs (|2.6|) have to be solved on the whole domain Q. Is it possible to 
localize the computation of the basis elements 9i to a neighborhood of the support of 
the elements fa? Observe that the support of each ci>; is contained in a ball B(x,,C h) of 
center X{ (the node of the coarse mesh associated with Xi) of radius C h. Let < a < 1. 
Solving the PDEs (|2.6p on sub-domains of O (containing the support of fa) may, a priori, 
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increase the error estimate in the right hand side of (|2.5|) , This increase can, in fact, be 
linked to the decay of the Green's function of the operator — div(oV). The slower the 
decay, the larger the degradation of those approximation error estimates. Inspired by 
the strategy used in [35] for controlling cell resonance errors in the computation of the 
effective conductivity of periodic or stochastic homogenization (see also [361 E31 [63] ) , we 
will replace the operator — div(aV) by the operator y — div(aV) in the left hand side 
of (|2.6p m order to artificially introduce an exponential decay in the Green's function. 
A fine tuning of T is required because although a decrease in T improves the decay of 
the Green function, it also deteriorates the accuracy of the transfer property. In order 
to limit this deterioration, we will transfer a vector space with a higher approximation 
order than the one associated with piecewise linear elements. Let us now give the main 
result. 

3.1 Localized bases functions. 

Let h € (0, 1). Let be an approximation sub- vector space of ifg(fi) such that 

• Xh is spanned by basis functions {tpi)i<i<N (with N = 0(\Q\/h )) with supports in 
B(xi,C h) where, the %{ are the nodes of a regular triangulation of O of resolution 
h. 

• Xh satisfies the following approximation properties: For all / £ Hq(Q) n H 2 (Q) 

JgJf-v\\ H t [n) <Ch\\f\\ HHn) , (3.1) 

and for all / G H%(Q) n H 3 (tt) 

Jg^ \\f-v\\ H i {a) < Ch 2 \\f\\ H3m . (3.2) 

• For all i, 

[ IV^I 2 < Ch d ~ 2 . (3.3) 
Jn 

• For all coefficients Cj, 

(3.4) 

i i 

Remark 3.1. Examples of such spaces can be found in |17j and constructed using piece- 
wise quadratic polynomials. From the first bullet point it follows that h can be though 
of as the diameter of the support of the elements ipi. The largest parameter h d /C satis- 
fying (|3.4p is the minimal eigenvalue of the stiffness matrix (J^(V</?j) T V^j)i<ij<Ar and 
Condition (j3.4j) is obtained from the regularity of the tessellation of 0. In fact, the 
proof of Proposition 13.21 shows that Condition (j3.4|) can be relaxed to the assumption of 
existence of a constant cL, > independent from h such that for all coefficients Cj 

(3.5) 

i i 
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Through this paper, we will write C any constant that does not depend on h (but 
may depend on d, J7, and the essential supremum and infimum of the maximum and 
minimum eigenvalues of a over fi). Let a G (0, 1) and C\ > 0. For each basis element 
<Pi of Xh let tpi be the solution of 




div(aVV>i) = Aifi in B(x h Cih a In \) n 9, 

on d(B(xi,Cih a ln^)nn). 



(3.6) 



Let 



14 := span(V>i) (3.7) 
be the linear space spanned by the elements ipi. 

Theorem 3.1. For g G L 2 (Q), let u be the solution of in Hq(£1) and the 

solution of (jl.ip in Vh- There exists Co > such that for C\ > Cq, we have 

-ufcllifi(n) ^ JCfc z/ aG (0,i] 



IMIl^) " [C/i 2 - 2a i/ a G [i,l), 
where the constants C and Cq depend on a, d, 6ut not on h. 

Remark 3.2. Theorem 13.11 shows the convergence rate in approximation error remains 
optimal (i.e., proportional to h) after localization if < a < 1/2 and decays to as 
/i 2 ~ 2q for i < a < 1. In particular, choosing localized domains with radii 0{^/h\nV) is 
sufficient to obtain the optimal convergence rate 0(h). Observe that the choice of the 
constant a in equation (|3.6|) is arbitrary. 

Remark 3.3. According to Theorem 13.11 the constant C\ in (|3.6p needs to be chosen 
larger than Co to achieve the convergence rate h + h 2 ~ 2a . The constant Co depends on 
a, d, A m i n (a) and A max (a). The constant C in the right hand side of (|3.8p also depends 
on a, d, A min (a) and A max (a). It is possible to give an explicit value for Co and C by 
tracking constants in the proof (in particular, as stated in Subsection l3.4l the dependence 
on A max (a) can be removed if the elements ^ are computed on sub-domains with added 
buffer zones around high-conductivity inclusions). 

Remark 3.4. If one uses piecewise linear basis elements instead of the elements (fi (i.e., 
in the absence of property (|3.2p ). then the estimate in the right hand side of (|3.8p 
deteriorates to /i 1-2a . The proof of this remark is similar to that of Theorem 13.11 The 
main modification lies in replacing h 2 /T by h/T in equations (|3.10p and (|3,16D . 

Remark 3.5. One could use piecewise linear basis elements instead of the elements <pi, 
and also remove the term h~ 2a ipi from the transfer property ()3.6|) . In this situation, 
we numerically observe a rate of convergence of h for periodic, stochastic and low- 
contrast media after localization of (|3.6p to balls of radii 0(h). In these particular 
situations (characterized by short range correlations in a), the term h~ 2a ipi should be 
avoided to obtain the optimal convergence rate h after localization to sub-domains of size 
0(h). In that sense, the estimate in the right hand side of (|3.8p corresponds to a worst 
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case scenario with respect to the medium a (characterized by long range correlations), 
requiring the introduction of the term and a localization to sub-domains of size 

0(\//iln^) for the optimal convergence rate h. 

Remark 3.6. For the elliptic problem, computational gains result from localization (the 
elements ipi are computed on sub-domains J7j of 0), parallelization (the elements ipi can 
be computed independently from each other), and the fact that the same basis can be 
used for different right hand sides g in (jl.ip . Computational gains are even more signif- 
icant for time-dependent problems because, once an accurate basis has been determined 
for the elliptic problem, the same basis can be used for the associated (parabolic and 
hyperbolic) time-dependent problems with the same accuracy (we refer to Sections H] 
and [5]). For the wave equation with rough bulk modulus and density coefficients, the 
proposed method (based on pre-computing basis elements as solutions of localized el- 
liptic PDEs) remains accurate, provided that high frequencies are not strongly excited 

(d t g e L 2 ). 

On Localization. We refer to [22], [25] and [6J for recent localization results for 
divergence- form elliptic PDEs. The strategy of [22] is to construct triangulations and 
finite element bases that are adapted to the shape of high conductivity inclusions via 
coefficient dependent boundary conditions for the subgrid problems (assuming a to be 
piecewise constant and the number of inclusions bounded). The strategy of [25] is to 
solve local eigenvalue problems, observing that only a few eigenvectors are sufficient to 
obtain a good pre-conditioner. Both [22] and [25] require specific assumptions on the 
morphology and number of inclusions. The idea of the strategy is to observe that if a 
is piecewise constant and the number of inclusions bounded, then u is locally H 2 away 
from the interfaces of the inclusions. The inclusions can then be taken care of by adapt- 
ing the mesh and the boundary values of localized problems or by observing that those 
inclusions will affect only a finite number of eigenvectors. 

The strategy of [BJ is to construct Generalized Finite Elements by partitioning the 
computational domain into to a collection of preselected subsets and compute optimal 
local bases (using the concept of n- widths [22] ) for the approximation of harmonic func- 
tions. Local bases are constructed by solving local eigenvalue problems (corresponding 
to computing eigenvectors of P*P where P is the restriction of a-harmonic functions 
from (jj* onto co C oj*, P* is the adjoint of P, and a; is a sub-domain of f2 surrounded by 
a larger sub-domain cj*). The method proposed in [6J achieves a near exponential con- 
vergence rate (in the number of pre-computed bases functions) for harmonic functions. 
Non-zero right hand sides (g) are then taken care of by solving (for each different g) par- 
ticular solutions on preselected subsets with a constant Neumann boundary condition 
(determined according to the consistency condition). 

As explained in Remark 12. 11 the near exponential rate of convergence observed in [B] 
is explained by the fact that the source space considered in [6] is more regular than L 2 
(since [UJ requires the computation particular (local) solutions for each right hand sides 
g and each non-zero boundary conditions, the basis obtained in [6] is in fact adapted to 
a-harmonic functions away from the boundary). The strategy proposed here can also 
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be used to achieve exponential convergence for analytic source terms g by employing 
higher-order basis functions <fi in (|3.6j) . Furthermore, as shown in sections HI l5l and 13.41 
the method proposed here allows for the numerical homogenization of time-dependent 
problems (because it does not require the computation of particular solutions for different 
source or boundary terms) and can be extended to high-contrast media. We also note 
that the basis functions tpi are simpler and cheaper to compute (equation (|3.6[) ) than 
the eigenvectors of P*P required by [6]. We refer to page 16 of [6] for a discussion on 
the cost of this added complexity. 

3.2 On Numerical Homogenization. 

By now, the field of numerical homogenization has become large enough that it is not 
possible to give an exhaustive review in this short paper. Therefore, we will restrict our 
attention to works directly related to our work. 

- The multi-scale finite element method IU2"l |4"T] can be seen as a numerical gener- 
alization of this idea of oscillating test functions found in .ff-convergence. A convergence 
analysis for periodic media revealed a resonance error introduced by the microscopic 
boundary condition [4Q|, 141] , An over-sampling technique was proposed to reduce the 
resonance error |40j . 

- Harmonic coordinates play an important role in various homogenization approaches, 
both theoretical and numerical. These coordinates were introduced in [42j in the con- 
text of random homogenization. Next, harmonic coordinates have been used in one- 
dimensional and quasi-one-dimensional divergence form elliptic problems [7J [5] , allowing 
for efficient finite dimensional approximations. The connection of these coordinates with 
classical homogenization is made explicit in [2] in the context of multi-scale finite ele- 
ment methods. The idea of using particular solutions in numerical homogenization to 
approximate the solution space of (jl.ip appears to have been first proposed in reservoir 
modeling in the 1980s [IB], [61] (in which a global scale-up method was introduced based 
on generic flow solutions i.e., flows calculated from generic boundary conditions). Its 
rigorous mathematical analysis was done only recently [19] and is based on the fact that 
solutions are in fact _ff 2 -regular with respect to harmonic coordinates (recall that they 
are i/ 1 -regular with respect to Euclidean coordinates). The main message here is that if 
the right hand side of (jl.ip is in L 2 , then solutions can be approximated at small scales 
(in ff 1 -norm) by linear combinations of d (linearly independent) particular solutions (d 
being the dimension of the space). In that sense, harmonic coordinates are only good 
candidates for being d linearly independent particular solutions. 

The idea of a global change of coordinates analogous to harmonic coordinates has 
been implemented numerically in order to up-scale porous media flows [27J ES HE] . We 
refer, in particular, to a recent review article [TBJ for an overview of some main challenges 
in reservoir modeling and a description of global scale-up strategies based on generic 
flows. 

- In |24] [29], the structure of the medium is numerically decomposed into a micro- 
scale and a macro-scale (meso-scale) and solutions of cell problems are computed on the 
micro-scale, providing local homogenized matrices that are transferred (up-scaled) to the 
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macro-scale grid. This procedure allows one to obtain rigorous homogenization results 
with controlled error estimates for non-periodic media of the form a(x, ~ ) (where a(x, y) 
is assumed to be smooth in x and periodic or ergodic with specific mixing properties 
in y). Moreover, it is shown that the numerical algorithms associated with HMM and 
MsFEM can be implemented for a class of coefficients that is much broader than a(x, -). 
We refer to |34] for convergence results on the Heterogeneous Multiscale Method in the 
framework of G and T-convergence. 

- More recent work includes an adaptive projection based method j3EJ, which is 
consistent with homogenization when there is scale separation, leading to adaptive al- 
gorithms for solving problems with no clear scale separation; fast and sparse chaos ap- 
proximations of elliptic problems with stochastic coefficients [60^ [371 [23] ; finite difference 
approximations of fully nonlinear, uniformly elliptic PDEs with Lipschitz continuous 
viscosity solutions [19] and operator splitting methods [H [3] . 

- We refer to |13[ [T2] (and references therein) for most recent results on homoge- 
nization of scalar divergence-form elliptic operators with stochastic coefficients. Here, 
the stochastic coefficients a(x/e,u) are obtained from stochastic deformations (using 
random diffeomorphisms) of the periodic and stationary ergodic setting. 



3.3 Proof of Theorem I3~TT! 

For each basis element ipt of X^, let ^i,r be the solution of 

1 7^i,T - div(aV^i,T) = A<£j in n 

The following Proposition will allow us to control the impact of the introduction of the 
term y in the transfer property. Observe that the domain of PDE (|3.9[) is still Q (our 
next step will be to localize it to f!j C fi). 

Proposition 3.1. For g £ L 2 {Q) let u be the solution of (jl.ip in Hq(Q,). Then, there 
exists v £ span('0j j r) such that 

\\u — v\\ H-ircn h 2 

] <C{h + %). (3.10) 



Furthermore, writing v := ^ (H^T we have 

Y,$<Ch- d {l + T- 2 )\\g\\l 2{n) (3.11) 

i 

Proof. Let v = ^ Cjf/^T- We have 

- div (aV(u - v ))=g+^-J2 d^Vi- (3-12) 



9 



Define a[v] to be the energy norm a[v] := J n (Vv) T aVv. Multiplying (|3TP21) by u — v and 
integrating by parts, we obtain that 



I l|2 

— — + a \u - v 



} = Jiu-vfa + ^-^aAtpi). (3.13) 



Write Cj = c^i + and let w\ and W2 be the solutions of Aw\ = g — J2i <k,iA<pi and 
Au>2 = If — Yli c i,2Afi with Dirichlet boundary conditions on d£l. Then, we obtain by 
integration by parts and the Cauchy-Schwartz inequality that 

II l|2 

II ^ ^11 1^- fn 1 ) 

+ a[u -v]< \\V{u - v)\\ {L2md (||Vwi|| (L 2 (n))(i + ||Vw2|| (i 2 (n))d ). (3.14) 

Using (|3.ip . we can choose (cf i) so that 

l|Vwi||( L 2 (n))d < Ch\\g\\ L 2 (n) . (3.15) 
Using (|3.2p . we can choose (c^) so that 

h 2 

ll v ^2||(L2 (c)) d < C—\\u\\ H i {n) , (3.16) 

we conclude the proof of the approximation (|3.10p by observing that < CIIs'IIl 2 ^)- 

Let us now prove Equation (|3.1ip . First, observe that Equation (j3.4j) and the triangular 
inequality imply that 

L 2 (n) + II ^Ci^ipiW^^Q)). (3.17) 



Next, we obtain from (|3.15|) and Poincare inequality and 

II E c *,i V ^Hi 2 (^) ^ C \\9\\l 2 {q.) (3-18) 

and 

II J2c i)2 V(pi\\ L 2 {n) < C-\\g\\ L2{n) (3.19) 

i 

We conclude by combining equations (|3.18p and (|3. 19|) with ()3.17p . □ 



We will now control the error induced by the localization of the elliptic problem (|3.9p . 
To this end, for each each basis element (fi of write Si the intersection of the support 
of (fi with fi and let fij be a subset of fi containing Si such that dist(5 , j, fi/fij) > 0. Let 
also V^T.fij be the solution of 

fy^TA - div(aVV'i,T,Oi) = Acpi in fi; 
[V'i.TA = on dfli. 

For i4, B C fi, write cZ(.A, B) the Euclidean distance between the sets A and -B. 
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Proposition 3.2. Extending ipi Tfti by on fi/fl$ we have 

ikt - *™n H , (n) < (di8t(S(i!i/S!j)r -> ( — c^—)- p' 21 ' 



We refer to Section [A] of the Appendix for the proof of Proposition 13.21 
Taking U { := B(x h Cxh a In ±) n (we use the particular notation C\ because our 
proof of accuracy requires that specific constant to be large enough, i.e., larger than 
a constant depending on the parameter C appearing in the right hand side of (|3.2ip 
and the parameter C describing the balls B(xi,C h) containing the support of the basis 
functions (</?i)i<i<jv introduced in Subsection 13. ip and T = h 2a in equation ()3.2ip of 
Proposition 13.21 we obtain for C\ large enough (but independent from h) that 

||^,T - ATMlHHn) ^ Ch d+1+2a . (3.22) 

Let u be the solution of (jl.ip in H^tt). Using Proposition 13. 1| we obtain that there 
exist coefficients Cj such that 

|h-^c^, T || H i (n) < C(h + h 2 - 2a )\\g\\ L 2 {n) . (3.23) 

i 

and 

^2 c "< Ch ~ d ~ 4a \\9\\Un) (3-24) 

i 

Using the triangle inequality, it follows that 
\\u - ^2 Ci^TA 11^1(0) - c { h + h2 ~ 2a )\\g\\LHn) + J }2\ci\\\^i,T - ^T.fii lljy^Q), (3.25) 

i i 
whence, from Cauchy- Schwartz inequality, 

\\u - ^tej^ip, <C(h + h 2 ~ 2a ) \\g\\ L 2 {n) 

i 

+ ( Yl l Cj l 2 ) 5 ( Yl IK T " feOillfrHO)) 2 • 
Combining ()3.26p with (|3,24p . we obtain that 
||« - Yl ^fen* llifi(a) + /l2 ~ 2 °) ll5lli2(Q) 



(3.26) 



+ C7i 2 2a ||s|| L 2(Q)(^||^,T-fenJ|#i 



(O) 



(3.27) 



Using (pT22]) in ([3T27]) . we obtain that 

Ih-E^feOiL^n) < c ( /l + /i2 " 2a )lbllL2 ( n ) . ( 3 _ 28) 

i 

Observe that it is the exponential decay in (|3.2ip that allows us to compensate for the 
large term on the right hand side of (|3.27|) via (|3.22p . This concludes the proof of 
Theorem 13.11 
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Figure 1: Illustrations of the buffer distance. 



3.4 On localization with high-contrast. 

The constant C in the approximation error estimate (|3.8p depends, a priori, on the 
contrast of a. Is it possible to localize the computation of bases for Vh when the contrast 
of a is high? The purpose of this subsection is to show that the answer is yes provided 
that there is a buffer zone between the boundaries of localization sub-domains and the 
supports of the elements ipi where the contrast of a remains bounded. More precisely, 
assume that Vt is the disjoint union of ^bounded and Qhigh- Assume that (|l,2p holds only 
on ^bounded, and that on VL high we have 



A mi „(a)|C| 2 <e T a(^<7lCP 



(3.29) 



where 7 can be arbitrarily large. Practical examples include media characterized by a 
bounded contrast background with high conductivity inclusions or channels. Let tp^ l9h 
be the solution of 



' h -2a^high _ div(aV ^) 
1>i = 



Aw 



in flj 
on 30, . 



(3.30) 



Let 

Vf 9h := span(Vf gh ) (3.31) 

be the linear space spanned by the elements ijj^ i9h . For each i, define b{ to be the largest 
number r such that there exists a subset 0,^ such that: the closure of contains the 
support of (fi, (£l'i) r is a subset of f2j (where A r are the set of points of that are at 
distance at most r for A), and (O-)V^ * s a subset of ^bounded- If n o such subset exists 
we set hi := 0. hi can be interpreted as the non-high-contrast buffer distance between 
the support of ipi and the boundary of Oj. We refer to Figure [TJ for illustrations of the 
buffer distance. 
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Theorem 3.2. For g G L 2 (Q), let u be the solution of in Hq{Q) and the 



where the constants C and Cq depend on A m i n (a), A max (a) (the bounds on a in ^bounded), 



Remark 3.7. Recall that the global basis computed in (j2.6|) remains accurate if the 
medium is both of high contrast (A max (a) >> 1) and degenerate (A m i n (a) << 1). The 
basis computed in (|3.30p preserves the former property (of accuracy for high contrast me- 
dia) but loses the latter (property of accuracy in the degenerate case) since the constant 
C in (|3.32p depends on A m i n (a). 

Remark 3.8. Observe that local solves have to resolve the connected components of high 
contrast structures. This is the price to pay for localization with high contrast in the 
most general case. Recall that in classical homogenization with high contrast the limit 
of the homogenized operator may be a non-local operator (we refer for instance to |21|). 
A similar phenomenon is observed here (distant points connected by high conductivity 
channels are associated with a low resistance metric and a large coupling coefficient in 
the numerically homogenized stiffness matrix). 

The proof of Theorem 13.21 is similar to that of Theorem 13.11 but it requires a precise 
tracking of the constants involved. Because of the close similarity we will not include 
the proof in this paper but only give its main lines. First, the proof of Proposition 
13.11 remains unchanged as the constants C in ()3.10|) and ()3.11|) do not depend on the 
maximum eigenvalue of the conductivity a. Only the proof of Proposition 13.21 has to 
be adapted and the part of the proof below Proposition 13.21 remains unchanged. This 
requires an application of the elements of lemmas IA.2} IA.31 IA.4I and IA.5I to buffer sub- 
domains (r^)V^i- The main point is to observe that the decay of the Green's function 
in (f2-) r /f^ can be bounded independently from 7 (due to the maximum principle). 

Observe that the choice of the sub-domain in (|3.30p can be chosen to be the same 
as in (|3.20p if its intersection with high contrast inclusions is void (i.e., if the maximum 
eigenvalue of a over Jlj remains bounded independently from 7); otherwise the choice 
of Qi in (|3.30p has to be enlarged (when compared to that associated with (|3.20p ) to 
contain the high-contrast inclusion (plus its buffer). 

4 The basis remains accurate for parabolic PDEs. 

The computational gain of the method proposed in this paper is particularly significant 
for time-dependent problems. One such problem is the parabolic equation associated 
with the operator — div(aV). More precisely, consider the time-dependent partial dif- 
ferential equation 




d, £1 but not on h and 7 (The upper bound on a on Qhigh)- 
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dtu(x,t) — div (a(x)Vu(x, t)\ = g(x,t) (x,t) G Sly; g G L 2 (Qt), 
u = on dQx, 



(4.1) 



where a and f2 satisfy the same assumptions as those associated with PDE (jl.ip . f2y := 
n x [0, T] for some T > and <9ft T := (<90 x [0, T]) U (fi x {t = 0}). 

Let V/j be the finite-dimensional approximation space defined in (|3.7p . Let Uh be the 
finite element solution of (|4.1|) . i.e., Uh can be decomposed as 

u h (x,t) =^ j c i (t)ij) i {x), (4.2) 

i 

and solves for all j 

(^j,d t u h ) L 2( n) = -a[ipj,u h ] + (ipj,g) L 2(n), (4.3) 
with a[v, w] : = J n (Vv) T aVw. Write 



NIlW,^)) := / \Vv\ 2 (x,t)dxdt. (4.4) 

»/ Jit 

Theorem 4.1. We have 

\\(u-u h )(.,T)\\ L2(n) + \\u - u h \\ L 2 (0 ^ H i {Q)) < C\\g\\ L 2 {nT) (h + h 2 - 2a ). (4.5) 



Proof. The proof is a generalization of the proof found in [50] (in which approximation 
spaces are constructed via harmonic coordinates). Let At be the bilinear form on 
L 2 (0,T,H^(n)) defined by 



At[wi,w 2 ]:= / a[wi,W2]dt. (4.6) 
J o 

Observe that for all v G L 2 (0, T, V h ), 

(v,d t {u - Uh)) L2 (Q T ) +A T [v,u- u h ) = 0. (4.7) 
Writing At[v] := •AtKv], we deduce that for u G L 2 (0,T, Vh), 

1 2 

-||(w-n h )(.,T)|| L2(f7) +^ T [u-^] = 
(it - u,<%(tt - ^h)) L 2 ( Q T ) + < At[«-'w,«- 
Using c^tt/j in (|4.3|) and integrating, we obtain that 

\\ d tUh\\h{n T ) + l a [uh(-,T),u h (.,T)} = {dtu h ,g) L2{UT y (4.9) 
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Using Minkowski's inequality, we deduce that 

\\dtu h \\ 2 L2{nT) + a[u h (.,T),u h (.,T)] < C|| 5 ||| 2(nr) . (4.10) 

Similarly, 

ll^ll£ 2( n T ) + a[u(.,T), u(, T)] < C\\g\\ 2 L2{nT y (4.11) 
Using Cauchy-Schwartz and Minkowski inequalities in (|4.8p , we obtain that 

ii 1 1 2 

|| (u - u fc )(.,T)|| L2(n) +A T [u - u h \ < C\\u - v\\ L 2 {nT) \\g\\ L 2^ T) + CA T [u - v). (4.12) 

Take v = TZhU to be the projection of u onto L 2 (0,T, V/J with respect to the bilinear 

form At- Observing that — div(aVu) = g — dtu with (g — dtu) € L 2 (J7r), we obtain 
from Theorem 13.11 that 

(At[u -K h u\) * < C\\g\\ L 2 {nT) (h + h 2 ~ 2a ). (4.13) 

Let us now show (using a standard duality argument) that 

||u - K h u\\ L 2 {nT) < C(h + h 2 - 2a ) 2 \\g\\ L 2 {nT) . (4.14) 

Choose v* to be the solution of the following linear problem: For all w £ L 2 (Q, T, Hq(Q)) 

At[w,v*] = (w,u-n h u) L 2(Q T y (4.15) 

Taking w = u — TZhU in (|4.15|) . we obtain that 

\\ u ~ Khu\\ 2 L 2( nT) = A T [u - TZ h u, v* - TZ h v*]. (4.16) 

Hence by Cauchy Schwartz inequality and (|4.13p . 

\\u - n h u\\ 2 L 2 (nT) < C{h + h 2 ~ 2a )\\g\\ L 2 (nT) (A T [v* - K h v*])K (4.17) 
Using Theorem 13.11 again, we obtain that 

(A T [v*-n h v*])* <c\\u-n h u\\ L 2 {nT) (h + h 2 ~ 2a ). (4.i8) 

Combining (gJHD with (|4"T?]) leads to (|4TT4"]) . Combining (|4TT2l with v = K h u, (|4TT4"]l 
and (|CT3D leads to 

|| (« - u h )(.,T)\\ 2 L2(n) + A T [u - u h ] < C(h + h 2 - 2a ) 2 \\g\\ 2 L2{QT) , (4.19) 

which concludes the proof of Theorem 14.11 □ 
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Discretization in time. Let (t n ) be a discretization of [0,T] with time-steps |i n +l ~~ 
t n \ = At. Write Zj., the subspace of L 2 (0,T, T4), such that 

Zj° = i.v € L 2 (0, T, V^) : = a(t)ipi(x), C{(t) are constants on (f n ,f n +i] j . 

(4.20) 

Write UhAu the solution in Zj. of the following system of implicit weak formulation 
(such that u/i i At( 2; 5 0) = 0): For each n and ip £ Vh, 

{ip, Uh,At(t n +i)) l2( q) ={i>i u h>At (t n )) L2(c) 

ftn + l (4.21) 

- \At\a[4),u hAt (t n+1 )] + (t/;, / g(t)dt) L2( py 

j t n 

Then, we have the following theorem 
Theorem 4.2. We have 



\{u - u hAt ){T)\\ L2{n) + \\u - u hAt \\ L 2 {QtTtH i m < C(\At\ + h + h 2 2a ) 

(ll%llL2(o,T,H- 1 (n)) + ||ff(-,0)|| L2(Q) + lbllL 2 (n T ) 



(4.22) 



The proof of Theorem 14.21 is similar to that of Theorem 1.6 of |50] and will not be 
given here. Observe that homogenization in space allows for a discretization in time 
with time steps 0(h + h 2 ~ 2a ) without compromising the accuracy of the method. 

5 The basis remains accurate for hyperbolic PDEs. 

Consider the hyperbolic partial differential equation 



p(x)d 2 u(x,t) - div \ a(x)Vu(x,t)j = g(x,t) (x,t) E g £ L (£1 T ), 

u = on dn T , ( 5 - 1 ) 
dtu = on ft x {t = 0}, 

where a, f2, Qt and 8£It are defined as in Section |U In particular, a is assumed to be 
only uniformly elliptic and bounded (ctjj € L°°(Q)). We will further assume that p is 
uniformly bounded from below and above (p 6 L°°(fi) and essinf p(x) > /? m i n > 0). It 
is straightforward to extend the results presented here to nonzero boundary conditions 
(provided that frequencies larger than 1/h remain weakly excited, because the waves 
equation preserves energy and homogenization schemes can not recover energies put 
into high frequencies, see [SI])- For the sake of conciseness, we will give those results 
with zero boundary conditions. PDE (|5.ip corresponds to acoustic wave equations in a 
medium with density p and bulk modulus a" 1 . 
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Let V/i be the finite-dimensional approximation space defined in (|3,7p . Let be the 
finite element solution of (|5.ip . i.e., Uh can be decomposed as 



Uh(x,t) = (x), (5.2) 

i 

and solves for all j 

{i>j,dtU h ) L 2 {Ptn) = -a[ipj,u h ] + (ipj,g)L*(n), (5-3) 

where 

(v>w)u*(p,n) / VW P- ( 5 - 4 ) 
Jn 

Theorem 5.1. If d t g G L 2 {Q, T ) and g(x,0) G L 2 {£1), then 
\\d t (u - u/j)(.,T)|| L2(n) + \\u - u h|| L 2 (0iTiH i( Q )) < 

C{\\d t g\\ L2{nT) + \\g(x,0)\\ L2{n) )(h + h 2 - 2a ). 



(5.5) 



Remark 5.1. We refer to |59j for an analysis of the sub-optimal rate of convergence 
associated with finite-difference simulation of wave propagation in discontinuous media 
(see also \18\ 156]). We refer to [51| for an alternative upscaling strategy based on har- 
monic coordinates. If the medium is locally ergodic with long range correlations [8] and 
also characterized by scale separation then we refer to HMM based methods \28\ Q]. 
Homogenization based methods require that frequencies larger than 1/h remain weakly 
excited. For high frequencies, and smooth media (or away from local resonances, e.g. 
local, nearly resonant cavities), we refer to the sweeping pre-conditioner method |3U|. I31| . 

Proof. Let At be the bilinear form on L 2 (0, T, Hq (fi)) defined in (|4.6|) . Observe that 
for all v G L 2 (0,T,V h ), 

(v,d?(u-u h )) L2[pnT) +A T [v,u-u h } = 0. (5.6) 

Taking dtu — dtUh — (dtu — dtv) as a test function in (|5.6p and integrating in time, we 
deduce that for d t v G L 2 (0,T,V h ), 

i||ft(«-ti h )(.,T)||J atan) + ia[(tt-« fc )(.,T)] = 
{ptiu - v), d 2 {u - u h )) L2{p QT) + A T [d t (u -v),u-u h ], 

where (v, w) L 2r p ^i T \ := J Q T JqVW pdxdt. Taking the derivative of the hyperbolic equa- 
tion for u in time, we obtain that 

d?u - div(aVd t u) = d t g. (5.8) 

Integrating (|5.8|) against the test function d 2 u and observing that d 2 u(x,0) = g(x,0), 
we also obtain that 

\\d 2 u(;T)\\l 2{pQ) +a[d t u(.,T)) < C(||%||| 2(Qt) + \\g(x,0)\\h { a))- (5-9) 
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h 


L 2 


H 1 


L°° 


0.5 


0.0119 


0.0913 


0.0157 


0.25 


0.0057 


0.0664 


0.0115 


0.125 


0.0027 


0.0482 


0.0075 


0.0625 


0.0005 


0.0207 


0.0032 



Table 1: Example 1 of Section 3 of [4jJ] (trigonometric multi-scale, see also [55]) with 
a = 1/2. 



Similarly, we obtain that 

\\d?u h (.,T)\\ 2 L2{p ^ + a[d t u h (.,T)] < C{\\d t g\\l HnT) + \\g(x,0)\\ 2 L2{n) ). (5.10) 

Take dtv = IZhdtu to be the projection of dtu onto L 2 (0,T, Vh) with respect to the 
bilinear form At- Observing that — div(aVdtu) = dtg — d 2 u with (g — d 2 u) £ L 2 (Q,t), 
we obtain from (|5.9p and Theorem 13.11 that 

{A T [u-n h u}f> <C(\\d t g\\ L 2 (QT) + \\g(x,0)\\ L 2 (n) )(h + h 2 ~ 2a ). (5.11) 

Furthermore, using the same duality argument as in the parabolic case, we obtain that 

||«-7^|| i2(pAr) < C(h + h 2 ~ 2a ) 2 {\\d t g\\ L 2 {nT) + \\g(x,0)\\ L 2 (Q))- (5.12) 

Using Cauchy-Schwartz and Minkowski inequalities and the above estimates in (|5.7[) , we 
obtain that 

\\d t (u-u h )(.,T)\\ 2 L2ipn) +a[(u-u h )(.,T)] < 

C(h + h 2 - 2a ){A T [u - u h ] + \\d t g\\ L 2 (nT) + \\g(x, 0)\\ L 2 (n) ) . 

We conclude using Grownwall's lemma. □ 



6 Numerical experiments. 
6.1 Elliptic equation. 

We compute the solutions of (jl.ip up to time 1 on the fine mesh and in the finite- 
dimensional approximation space Vh defined in (|3.7p . The physical domain is the square 
[— 1, l] 2 . Global equations are solved on a fine triangulation with 66049 nodes and 131072 
triangles. 

The elements {ipi) of Sub-section 13. II are weighted extended B-splines (WEB) [381139] 
(obtained by tensorizing one-dimensional elements and using weight function (1— j; 2 )(1 — 
y 2 ) to enforce the Dirichlet boundary condition). The order of accuracy is not affected 
by the choice of weight function given that the boundary is piecewise smooth. Our 
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1 1.5 2 

(a) L 2 error 



-> -5 L 
2.5 0.5 



1.5 2 

(b) H 1 error 



2.5 



Figure 2: Example 5 of Section 3 of |49j (percolation at criticality) . Logarithm (in base 
2) of the error with respect to log 2 (/io//i) (for h = 0.125) and the value of T used in 



motivation for using WEB elements lies in the fact that, with those elements, (Dirichlet) 
boundary conditions become simple to enforce. This being said, any finite elements 
satisfying the properties (|3.ip . (|3.2|) . (|3.3p and (|3.4j) would be adequate [II] . 

We write h the size of the coarse mesh. Elements ipi are obtained by solving (|3.6|) 
on localized sub-domains of size ho- Table [I] shows errors with a = 1/2 for a given by 
(|6.ip (Example 1 of Section 3 of [49] . trigonometric multi-scale, see also [45]). i.e., for 



1 l.l + sin(2vrx/ei) 1.1 + sin(2vry/e 2 ) 1.1 + cos(2vrx/e 3 ) 

a ^ X ' + sin(2vry/ei) + 1.1 + cos(2vnr/e 2 ) + 1.1 + sin(27n//e 3 ) + 

1.1 + sin(2vry/e 4 ) 1.1 + cos(27rx/e 5 ) , 2 2un 
H r-7 / — r + sin(4a; y ) + 1), 



(6.1) 



1.1 + cos(2-7rx/e4) 1.1 + sin(27ry/e5) 

where ei = ± e 2 = ^,e 3 = ^,e 4 = ^,e 5 = ^. 

Figure [2] shows the logarithm (in base 2) of the error with respect to log 2 (/io//i) (for 
h = 0.125) and the value of T used in (|3.6p for a given by Example 5 of Section 3 of 
[4"9"] (percolation at criticality, the conductivity of each site is equal to 7 or I/7 with 
probability 1/2 and 7 = 4). 

Figure [3] shows the logarithm (in base 2) of the error with respect to log 2 (/io//i) (for 
h = 0.125) and the value of T used in (|3.6p for a given by Example 3 of Section 3 of [49j, 
i.e., a{x) = e h ( x \ with h(x) = X^|fc|<ij( a fc sin(2-7r/c • x) + cos(27r/c • x)), where and bk 
are independent uniformly distributed random variables on [—0.3, 0.3] and R = 6. 

Remark 6.1. Two factors contribute to the error plots shown in figures [2] and [3] a 
localization error which becomes dominant when ho/h is small (i.e. the fact (|2.6|) is not 
solved over the whole domain Q) and the distortion of the transfer property resulting 
from the 1/T term in (|3.9p . As expected both figures show that when ho/h is large, 
the error due to the distortion of the transfer property is dominant and is minimized 
by a large T. However, when ho/h is small, the localization error is dominant and 
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1 1.5 2 

(a) L 2 error 



-J _4L 
2.5 0.5 



1.5 2 

(b) H 1 error 



2.5 



Figure 3: Example 3 of Section 3 of [49] (exponential of a sum of trigonometric functions 
with strongly overlapping frequencies). Logarithm (in base 2) of the error with respect 
to log 2 (/io/M (f° r h = 0.125) and the value of T used in ()3.6j) . 




Figure 4: High conductivity channel. 



is minimized by a small T. The fact that in Figure [3] this error is minimized by the 
second smallest T instead of the smallest T is explained by the fact that the localization 
error remains bounded when ho/h is of the order of one whereas the error due to the 
distortion of the transfer property blows up as T ! 0. The fact that in both figures, 
curves associated with different T interest each other, is indicative of the fact that 
for intermediate values of ho/h, the error can be minimized via a fine-tuning of T as 
explained in section The differences in the locations of these intersections can be 
explained by a larger localization error associated with the example of Figure [3] (due 
to longer correlation ranges). In particular, the comparison between figures [5] and [3] 
indicates larger errors for Figure 
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(a) L 2 error (b) H 1 error 



Figure 5: High conductivity channel (Figure U|). The x-axis shows log 2 (/i), the y-axis 
shows the log 2 of the error in I? and i/ 1 -norm. The three cases for the localization are 
ho = 0(y/h\n^) with a buffer around the high conductivity channel (see Sub-section 
I3.4p of size 0(\/7im ho = 3h with no buffer around the high conductivity channel 
and ho = 3h with a buffer around the high conductivity channel of size 3h. 



6.2 High contrast, with and without buffer. 

In this example, a is characterized by a fine and long-ranged high conductivity channel 
(Figured]). We choose a(x) = 100, if x is in the channel, and a(x) is the percolation 
medium, if x is not in the channel (the conductivity of each site, not in channel, is equal 
to 7 or I/7 with probability 1/2 and 7 = 4). Figure [5] shows the log% of the numerical 
error (in L 2 and H 1 norm) versus log2(h). The three cases for the localization are 
ho = 0(^/hln j-) with a buffer 6, around the high conductivity channel (see Sub-section 
13, 4p of size 0{-\fhhi h), ho = 3/i with no buffer around the high conductivity channel and 
ho = 3/i with a buffer 6j around the high conductivity channel of size 3h. The first case 
shows that the method of Sub-section 13.41 is converging as expected. The second case 
shows that, as expected, taking a = 1, does not guarantee convergence. The third case 
shows that adding a buffer around the high conductivity channel improves numerical 
errors but is not sufficient to guarantee convergence (as expected, we also need a < 1). 
The percolating background medium has been re-sampled for each case; the effect of this 
re-sampling can be seen for the largest value of h (i.e., log 2 (/i) = —1). 

6.3 Wave equation. 

We compute the solutions of (|5.ip up to time 1 on the fine mesh and in the finite- 
dimensional approximation space Vh defined in (|3.7p . The initial condition is u{x, 0) = 
and ut(x,0) = 0. The boundary condition is u(x,t) = 0, for x £ d£l. The density is 
uniformly equal to one and we choose g = sin(-7ra;) sin(7ry). Figure [6] shows the fine mesh 
solutions u and at time one, for a given by the trigonometric example (|6.1|) . with 
h = 0.125, ho = 3h and T = h. Figure [6] shows the fine mesh solutions u and Uh at time 
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(a) u (b) Uh 

Figure 6: Wave equation. Trigonometric case, fine mesh solution, h = 0.125, ho = 3h, 
T = h. The L 2 , H 1 and L°° relative numerical errors are 0.0339, 0.1760 and 0.0235. 



one, for a given by the high conductivity channel example (Figure U]), with h = 0.125, 
ho = 3h and T = h. 

We refer to |52] for a list of movies on the numerical homogenization of the wave 
equation with and without high contrast and with and without buffers (extended buffers 
in the high contrast case). 



A Proof of Proposition 13.2 



The proof of Proposition l3.2l is a generalization of the proof of the control of the resonance 
error in periodic medium given in [35]. 

First we need the following lemma, which is the cornerstone of Cacciopoli's inequality. 



(A.l) 



Lemma A.l. Let D be a sub-domain of T2 with piecewise Lipschitz boundary, and let v 
solve 

j f, - div (a(x)Vv(xfj = f(x) xeD;f e H~ 1 (D), 
\v = on dD, 

Let £ : D — > M + be a function of class C 1 such that C is identically null on an open 
neighborhood of the support of f . Then, 

[ \V(Cv)\ 2 <C [ v 2 \VC\ 2 , (A.2) 

J D JD 

where C only depends on the essential supremum and infimum of the maximum and 
minimum eigenvalues of a over D. 

Proof. Multiplying (jA.ip by C^v and integrating by parts, we obtain that 

J + j D V{t 2 v)aVv = Q. (A.3) 
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(a) u (b) u h 

Figure 7: Wave equation. Channel case, coarse mesh solution, h = 0.125, ho = 3h, 
T = h. The L 2 , H 1 and L°° relative numerical errors are 0.0439, 0.2684 and 0.0389. 



Hence, 

/ C V - + I V((v)aV((v) = ( v 2 V(aV(, (A.4) 

JD 1 JD JD 

which concludes the proof. □ 

Lemma A. 2. Let D be a sub-domain of Q with piecewise Lipschitz boundary. Write 
Gt,d the Green's function of the operator y — div(aV) with Dirichlet boundary condition 
on dD. Then, 

GtM^v) < y^z^pi ex p(-^r7=b> ( A - 5 ) 



where C only depends on d and the essential supremum and infimum of the maximum 
and minimum eigenvalues of a over D. 

Proof. Extending a to M d and using the maximum principle, we obtain that 

Gt,d(x,v) < G T>R d(x,y), (A. 6) 

we conclude by using the exponential decay of the Green's function in R d (we refer to 
Lemma 2 of E5]). □ 



Lemma A. 3. Let ipi? be the solution of (|3.9p and the solution of (|3.20p . Let 

Q'i be a sub-domain o/Oj such that Si C and dist(5j, Cli/Q'j) > 0. We have 

11 , 11 Chi~ l f dist(5i,fi/0')\ 



and 



1 , I, chi~ l f dist(5i,n/n;-) . 

"«* "'"«"■'"■» - (i*t(s„n/n d y exp y — c&— 1 • 
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Proof. For A C f2, write A r the set of points of O that are at distance at most r from 
A. Let us now use Cacciopoli's inequality to bound Jq/q/ |VV>j,T| 2 - Using Lemma [A. II 
with C identically equal to one on f2/^, zero on (0/^) r with r := dist(Sj, and 
|VC| < C/r, we obtain that 



/ 

Jn, 



m/n 1 . r - J(si/sitf 

Next, observe that for x £ (fi/0^) r 



|V^,t| 2 < § / $, T . (A.9) 



A,t( x ) = ~ VG T ,ti{x,y)V<Pi(y)dy. (A.IO) 
■/si 

Hence, 

|^i,r(aO| < ||Vy>i||( L 2(5 i ))d||VGr,n(x,.)|| (La(S4)) d- ( A - n ) 
Another use of Cacciopoli's inequality leads to 

C 

HVGt.q^, 011(^2(5.)). < -\\Gt,q{x, .)\\ L 2 {S ry (A.12) 
Combining (|A9j) with (fATTj) with (jA~T2"|) . we obtain that 

L K m >'* 2 "^"kw l wa , r II^C ->ll W (A.i3) 

We conclude the proof of ()A.7p using Lemma [A.2l and (|3.3p . The proof of (|A.8j) is similar 

observing that dist(S , i , Q/ft-) < dist(S^, fit/ft-) □ 

Lemma A. 4. Let D be a sub-domain of £1 with piecewise Lipschitz boundary. Let 
ip £ i7 1 (f]), and Zei v sofee 



'f, -div (o(x)V!)(i)j 
v = tp on dD, 



x G L>, , , , 

(A.14) 



Write S the intersection of the support of ip with D. Let D\ be a sub-domain of D such 
that dist(Z?i,5') > 0, then 

L |Vt|2 5 i^±w AT ^ + lfm ™ exp ( - (A - 15) 

where C does not depend on D,D\,S. 
Proof. Write w := v — tp. Then, 



'f - div (a(x)Vw(x)) = -% + div(aVV>) x £ D, 
v = on 3D, 



(A.16) 
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Thus, 

w(x) = - J (^G TiD (x, y) + V^(y)a(y)VG T:D (x, y)) dy. (A.17) 
Using Cauchy- Schwartz inequality, we obtain that 

\w(x)\ < C||V||^i(c)(^||G T , D (x,.)|| L 2 (s) + ||VG T ,D^,.)|| (L 2 {s))d )- (A.18) 

For A C D, write A r the set of points of D that are at distance at most r from A. 
Let us now use Cacciopoli's inequality to bound J D |Vu;| 2 . Using Lemma lA.ll with £ 
identically equal to one on D\, zero on D/D 7 ^ 1 and such that |V£| < Cjr\ we obtain 
that 

/ |VH 2 <4 [w 2 , (A.19) 
Jd 1 r i JdI 1 

provided that dist(D[ 1 , S) > 0. Hence, for r\ := distfiJi . S)/3, we obtain ()A,19p . Taking 
T2 := dist(Z>i, S)/3 and using Cacciopoli's inequality again, we also obtain that 

C 

\\VG t ,d{x, .)|| (L 2 (5))d < — \\Gt,d{x, -)|L2 (S r 2) . (A.20) 

Combining ()A.19p with ()A.18j) and (|A.20j) and observing that w = v on D^ 1 we obtain 
that 

JDi 'VI JD^ 

Using Lemma lA.21 we deduce that 

k ' v '"' 2 * (d^s)H *'"<" (rl + I,2exp( " dJ ^ y (A ' 22) 

This concludes the proof of Lemma IA.41 □ 

Lemma A. 5. Let ipi? be the solution of (|3.9|) and ipi t T,tti the solution of (|3.20|) . Let 
be a sub-domain of fij such that dist(f2/fij, > 0. VFe /icwe 

1 1 , „ CfT- 1 + - 1 ( dist(^M,0')\ /A s 

Proof. Lemma IA.5I is a direct consequence of Lemma IA.41 To this end, we choose 
D := Qi, v =: ij)i T~ "0iTQj and D\ := f2^. We also choose ip := tji/j^t where rj : f2 — )■ [0, 1] 
is C 1 , equal to one on fi/f2j and on (f2/T2j) r with r := dist(f2/fij, f^)/3 (j4 r being the 
set of points in Q at distance at most r from A) and |Vr/| < C/r. We obtain from 
Lemma IA.41 that 



;,T 



( J (dist (n/n^n ) V CVT J 



We conclude using and ||^||iri(n) < Sitpfen^ ll ^11(^(0))"- D 
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Observing that 



+ \\ll>i,T 



(A.25) 



we conclude the proof of Proposition ^. 2l bv using Lemma [A.5l and Lemma [A.3l with := 
S\ where S\ are the points in fij at distance at most r from Si with r := dist(5j, Q/Cli)/3. 

B On the compactness of the solution space. 

Although the foundations of classical homogenization [U] were laid down based on as- 
sumptions of periodicity (or ergodicity) and scale separation, numerical homogenization, 
as described here, is independent from these concepts and solely relies on the strong com- 
pactness of the solution space (and the fact that a compact set can be covered with a 
finite number of balls of arbitrary sizes). Observe that an analogous notion of com- 
pactness supports the foundations of G and .ff-convergence ( [47] . \57\ I32j). The main 
difference is that G and iif-convergence rely on pre-compactness and weak convergence 
of fluxes and here, we rely on compactness in the (strong) f^Q-norm, i.e. the following 
theorem. 

Let W be the range of g in . Write 



Theorem B.l. Let v < 1. If W is a closed bounded subset of H "(O) then W is a 
compact subset of Hq(Q) (in the strong Hy-norm). 

Proof. We have (aVu) po t = -VA _1 g. So using the same notation as in (|2.4|) we get 
(aW)p t = — VA" 1 ^. Let u n be a sequence in V then there exists a sequence in 
W such that — div(aVit n ) = g n . Using the fact that —V/S.~ l W is a compact subset 
of (L 2 (£l)) d (we refer, for instance, to the Kondrachov embedding theorem) we get 
that there exists g* G W such that ||VA _1 5f n - X7A~ 1 g*\\ L 2 — > 0. Writing u* the 
solution of — div(aVti*) = g* and using (aV(n n — u*)) po t = —VA~ 1 (g n — g*) we get that 
|| (aV(n n — u*)) po t\\L2 — > 0. Using the equivalence between the flux norm and the Hq 
norm we deduce that \\u n — u*\\ H i — > 0. This finishes the proof. □ 

This notion of compactness of the solution space constitutes a simple but fundamental 
link between classical homogenization, numerical homogenization and reduced order 
modeling (or reduced basis modeling [221 S3]) (we also refer to the discussion in Section 
6 of |10|). This notion is also what allows for atomistic to continuum up-scaling [63], the 
basic idea is that if source (force) terms are integrable enough (for instance in 1? instead 
of H~ l ) then the solution space is no longer H 1 but a sub-space V that is compactly 
embedded into H 1 and, hence, it can be approximated by a finite-dimensional space (in 
i? 1 -norm). In other words if these systems are "excited" by "regular" forces or source 
terms (think compact, low dimensional) then the solution space can be approximated by 
a low dimensional space (of the whole space) and the name of the game becomes "how to 



V := {u G Hq{VL) : u solves (fTTj) for some g G W}. 



(B.l) 
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approximate" this solution space (and this can be done by using local time-independent 
solutions). 
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